AD-A064  266 


UNCLASSIFIED 


UNIVERSITY  OF  SOUTH  FLORIDA  TAMPA  DEPT  OF  ELECTRICAL— ETC  F/G  13/10,1 
DESIGN  OF  FEEDBACK  CONTROL  AND  GEOMETRY  PARAMETFRS  VIA  MOFNM.(U) 

OCT  78  V K JAIN*  G L LOGSTON  N61339-75-C-0122 


SS-206-FB 


NL 


I » I 

ADA 

064266 


END  END 


DATE  DATE 

FILMED  FILMED 

4 -79  4 -79 

DDC  DOC 


Engineering  Research  Report 
No.  SS-20G+FB 


& 

CO 

CM 

CO 


t 


Q_ 

CD 

CD 


DESIGN  OF  FEEDBACK  CONTROL 
AND  GEOMETRY  PARAMETERS  VIA  flOFHM 

V.  K,  JAIN 
G.  L.  LOGSTON 

Department  of  Electrical  Engineering 
University  of  South  Florida 
Tampa,  Florida  33620 


OCTOBER  1978 

Approved  for  public  release;  distribution  unlimited 


Prepared  for: 

NAVAL  COASTAL  SYSTEMS  CENTER 
PAN  AM, A CITY,  FLORIDA  32RC1 


Contract  No.  N61 339-75-C-01 22-HR38 


DISCLAIMER  NOTICE 


THIS  DOCUMENT  IS  BEST  QUALITY 
PRACTICABLE.  THE  COPY  FURNISHED 
TO  DDC  CONTAINED  A SIGNIFICANT 
NUMBER  OF  PAGES  WHICH  DO  NOT 
REPRODUCE  LEGIBLY. 


1 


Unclassified 


N OF  THIS  PAGE  W idn  Data  Entered) 


READ  INSTRUCTIONS 
BEFORE  COMPl  ETINC  FORM 


REPORT  DOCUMENTATION  PAGE 


1 REPOHl  NUMBER 


2.  GOVT 


SS-20G+FB 


S.  TYPE  OF  REPORT  4 PERIOD  COVEREO 


4.  TITLE  (and  Subtitle) 


Design  of  Feedback  Control  and  Geometry 
Parameters  via  MOFNM  - 


C.  PERFORMING  GAG.  REPORT  NUMBER 


7.  AUTHOR^*; 


N6.1339-75-C 


Jain 


9.  PERFORMING  ORGANIZATION  NAME  ANO  AOORESS 


Department  of  Electrical  Engineering 
University  of  South  Florida 
Tampa,  Florida  33620 


P.E.63561N,  Pffly207 ■ 

T. A. 50207,  kUO. 20812-1 


11.  CONTROLLING  OFFICE  NAME  AND  ADDRESS 

Naval  Coastal  Systems  Center 
Panama  City,  Florida  32407 


September  1978 


is.  SECURITY  C_»SS.  (ol  tnta  teport) 


14.  MONITORING  AGENCY  NAME  A ADDRESSm  dllleient  Irom  Central 


Unclassif led 


16  DISTRIBUTION  ST  A"  £mEN  T (ol  IhU  Report) 


Approved  for  Public  Release;  Distribution  Unlimited 


17.  DIST  HIB'JTICN  STATEMENT  (ol  Ihr/cbaUa^u.  tn  BlbcJ.  It  <J/71er»rir  fretm  Haport) 


18  SUPPU'MEnTARY  NOTES 


19  KEY  WORDS  'Contlnu*  on  rovara*  aide  It  nacaaaary  and  Idantlly  by  block  numbar) 


Adaptive  Optimization 

Optimal  Design 

Desired  Trajectory 

Design  Trajectory  Error  Function 


Geometry 

Feedback  Parameter  Optimization 
Multiple  Object  Function-Newton  Method 
Submerged  Vehicle  Control 


i ra*ataa  aide  II  nacaaaary  and  Idantlly  by  block  numbar) 


ABSTRACT  (Contlnua  a 

^ A modified  Newton  technique  (MOFNM)  is  presented  and  applied  to  the 
optimization  of  several  design  parameters  of  a remotely  piloted  vehicle. 
The  technique  is  programmed  to  simultaneously  optimize  the  gains  of  the 
feedback  system  and  the  geometry  parameters  of  the  vehicle.  A logrithmic 
transformation  is  developed,  and  implemented  to  guarantee  sign  definite 
final  design  parameters.  Several  examples  are  included  to  illustrate  the 
technique's  effectiveness,  -yr  / / /z  s-  — 


EDITION  OF  ' NOV  CS  IS  OBSOLETE 
S N 0 102-  LF-  014-  66P1 


Unclassi f led 

security  classTfication'of  Thi 


Table  of  Contents 


Section  Title  Page  No. 


I Theory  1 

II  Feedback  Parameter  Optimization  4 

A.  Optimal  Design  with  Pitch  and  Depth  Command  Doublet  Input  4 

B.  Optimal  Design  with  Pitch  Command  Step  Input 

III  Geometry  and  Feedback  Optimization  (Pitch  Maneuvers)  7 

A.  System  Design  with  Factor  =0  * 7 

B.  Design  using  Adaptive  Method  9 

C.  Design  using  J0PT2  Adjustment  14 

IV  Feedback  and  Geometry  Optimization  (Depth  Maneuvers)  18 

A.  Optimal  Design  with  Ramp-Step  Input  19 

B.  Optimal  Design  with  Triplet  Input  26 

V Parameter  Optimization  Utilizing  Exponential  Transportation  33 

Conclusions  38 

References  39 

Appendix  A 40 

Appendix  B 42 

Appendix  C 43 

Appendix  D 47 

Appendix  E 49 


List  of  Figures 


Figure  No. 


Page  No. 


Generalized  Control  System 
Stern  Plane  Feedback  System 


2 

2 


3 

Comparison  of  Desired 
Doublet  Input 

and 

Final 

Responses , 

6 

4 

Comparison  of 
Step  Input 

Desired 

and 

Final 

Responses , 

8 

5 

Comparison  of 

Desired 

and 

Final 

Responses,  P=0 

10 

6 

Parer  vs  Factor 

11 

7 

Parer  vs  P 

11 

8 

Comparison  of  Desired  and 
Example  B-l  Adaptive 

Final 

Responses , 

13 

9 

Comparison  of 
Example  B-2 

Desired  and 
Adaptive 

Final 

Responses, 

15 

10 

Comparison  of  Desired 
J0PT2  Adjustment 

and 

Final 

Responses , 

17 

11 

12 

Depth  Command 

Comparison  of 
Example  A-l 

(z  ) 

com 

Desired 

and 

Final 

Responses , 

19 

21 

13 

Comparison  of 
Example  A-2 

Desired 

and 

Final 

Responses , 

23 

14 

Comparison  of 
Example  A- 3 

Desired 

and 

Final 

Responses , 

25 

16 

Depth  Command 

Input  (Z  ) 

com 

26 

17 

Comparison  of 
Examp le  B- 1 

Desired 

and 

Final 

Responses, 

28 

18 

Comparison  of 
Example  B-2 

Desired 

and 

Final 

Responses , 

30 

19 

Comparison  of 
Example  B-3 

Desired 

and 

Final 

Responses , 

32 

20 

Comparison  of 
Example  1 

Desired 

and 

Final  Responses, 

35 

21 

Comparison  of 
Example  2 

Desired 

and 

Final  Responses, 

37 

r 


DESION  OK  FEEDBACK  CONTROL  AN.)  CEOMETRY  PARAMETERS  VIA  MOKNM 


A significant  breakthrough  in  geometry  design  and  optimization  was  achieved 
in  the  research  effort  of  1976-77  [ 1 ] . Briefly,  a modified  Newton  technique 
(Multiple  Object-Function  Newton  Method)  was  developed  and  programmed  to  perform 
simultaneous  optimization  of  several  geometry  parameters,  e.g.,  control  surface 
areas,  lengths,  etc.  The  potential  benefits  include  improved  deflector  shapes 
for  existing  hull  designs,  and  reduction  of  design-evaluation-redesign  cycle 
time  for  completely  new  crafts.  Pursuing  this  work,  two  further  advances  are 
attempted  in  the  present  work.  First,  gains  of  the  feedback  control  system  and 
the  geometry  parameters  are  collectively  considered  for  optimization.  The 
success  of  this  effort  should  bring  about  closer  collaboration  between  the 
body-shape  designer  and  the  control  system  engineer.  Even  configurations  not 
considered  heretofore  could  be  evaluated  rapidly  for  their  true  performance 
potential  — and  used  when  deemed  superior  by  the  engineering  team.  A second 
improvement  considered  is  in  the  mathematical  formulation  of  the  optimization 
problem.  By  use  of  a logarithmic  transformation,  the  resulting  computer  solution 
is  sought  to  be  sign  definite,  and  is  thereby  guaranteed  to  be  physically  realizable. 

In  summary,  a methodology  for  the  designer  is  now  available  so  that  tie  may 
harness  the  full  potential  of  body-shape  — including  deflection  surfaces  — and 
control  gains  for  maximum  performance  of  the  craft. 

Examples  presented  here  pertain  only  to  the  longitudinal  dynamics. 


1.  THEORY 

The  linearized  state  equations  of  a vehicle*  are  of  the  form  [3],  [5] 

A — x = Bx  + Cu  (1) 

dt 

We  will  assume  that  the  longitudinal  and  lateral  dynamics  can  be  considered 
decoupled,  [4]  and  thus  can  be  analyzed  independently.  Concentrating  then  on 
longitudinal  dynamics,  equation  (1)  can  be  used  to  characterize  the  response  of 
variables.  Specifically,  the  state  vectors  become 

(2) 


* The  vehicle  under  consideration  is  a remotely-piloted  vehicle  (RPV) : It's 
hydronamic  coefficients  for  longitudinal  dynamics  are  listed  in  Appendix  I). 
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while  the  control  deflection  vector  is 


The  control  system  configuration  shown  in  Figure  1 results  in  the  feedback  law 


Figure  1.  Generalized  Control  System 


In  this  phase  of  the  study,  we  will  use  the  stern  plane  as  the  only 


control  input.  The  feedback  configuration  used  is  shown  in  figure  2. 


See  Apendix  A for  the  definition  of  other  design  parameters. 


2 


The  feedback  law  which  governs  the  system  depicted  in  figure  2 is  given  by 


u = 6 


-K®  Ks  COm  + 

° ZJ  Z 

com  _ 


0 - K?  Ka  (K„K®  + K?U  )-  K KS 
z o 2 0 z o lz 


Given  the  design  responses  z(k) ,k=l,2> . . . ,K,  estimates  for  the  optimum 
feedback  and  geometry  parameters  R are  found  such  that 

J = l (z(kA)-x(kA)  ]TQ[z(kA)-x(kA)  ] + [R0-RvfllTP 
k=l 

is  minimized  [6],  [8]  where 

A A A A A rp 

R = [ K.,K-  K K X ] 
z U z U s 

Equation  (6)  can  be  written  [1]  as 
£ 

J =J:X  tz(kA)-xv(kA)-H1H2(%n“V]TQ[z(kA)"xv(kA)‘HiH2(Vr%] 


u 

W 

6 (5) 

0 

Z 


where 


+ [VVii  p iRo-Vi] 


H,  = -5* 

1 3CT 

H = ~ 

2 3RT 

c = f(R)  Hydrodynamic  coefficients 

A 

Setting  the  partial  derivative  with  respect  to  R equal  to  zero,  [2],  [12]  we 


ob  t a in 


|j[  = 0 = -2  l H2rH]TQ[z(kA)-xv(kA)-H1H2(Rv+1-Rv)] 
k=l 

-2RTPR(Ro-Rv+1) 

where  R is  a diagonal  scale  matrix  for  the  feedback  and  geometry  parameters 

A 

R such  that 


R = RR 


A 

The  solution  to  (12)  for  the  new  value  of  the  feedback  parameters  R is  given  by  [1] 


(13) 


The  optimization  method  oulined  above  is  used  in  the  multiple  object  function 
approach,  MOFNP,  along  with  the  constraint  (banded  prediction) 

x(kA  + A)  = c(k)[z(kA)  + AA_1  (B  z(kA)  + Cu(kA))] 

+ (l-C(k))  [x(kA)  + AA-1  (B  x(kA)  + Cu(kA))]  (u) 

x(0)  = x°  initial  conditions 

where  £(k)  is  an  appropriately  chosen  sequence  of  0's  and  l's. 

II.  FEEDBACK  PARAMETER  OPTIMIZAION 

This  section  deals  with  the  selection  of  the  feedback  parameters  of  figure  2 
that  will  yield  the  optimum  trajectory  with  respect  to  a specified  desired 
trajectory.  The  stern  plane  geometry  parameter  A^  has  been  hard-wired  to  1.0,  thus 
eliminating  geometry  parameter  optimization  for  the  present. 

A . Optimal  Design  with  Doublet  Input 

Consider  the  system  given  in  figure  2,  excited  with  the  foliowing  input 
combination : 

i)  1 degree  pitch  angle  doublet  with  12.5  seconds  positive  and  12.5 
seconds  negative. 

ii)  100  feet  depth  command  doublet  with  12.5  seconds  positive  and  12.5 
seconds  negative. 

It  can  be  shown  that,  in  order  to  assure  system  stability,  the  following 
conditions  must  be  met:  1)  The  feedback  gains  R^  ^ and  R^ , (see  figure  A1  in 
Appendix  A)  must  provide  unity  feedback  in  order  to  generate  the  actual  pitch 
and  depth  error  signals.  2)  Since  the  numerical  value  of  the  depth  command 
is  two  to  three  orders  of  magnitude  larger  than  the  pitch  angle  command.  Kg 
should  be  two  to  three  orders  of  magnitude  larger  than  in  order  to  assure 

that  comparable  contributions  to  the  control  input  are  produced.  3)  The  depth 
rate  feedback  gain  K^  must  be  chosen  to  be  small  to  prevent  excessive  overshoot 
and  ringing  in  the  response.  4)  The  pitch  rate  feedback  gain  K/,  was  chosen  to 
be  about  two  orders  of  magnitude  greater  than  to  control  the  pitch  rate. 

These  conditions  were  applied  to  the  optimization  of  the  stern  plane  feedback 
parameters , 


- 4 


t 


R ■ [v*ii'V*< 


where  A , the  stern  plane  geometry  parameter,  was  hard-wired  to  1.0.  The 
results  of  the  experiment  are  given  in  table  1;  the  program  settings  are  given 
in  Table  2. 

Table  1.  Feedback  Parameters  and  Errors  , Doublet  Input  Response 


Feedback  Parameters 

VV  S'V  W K0<R12! 

RMS  % Difference  (%) 

U W 0 0 Z 

Baseline 

0.02 

1.80 

0.05 

9.50 

75.8 

mm 

67.9 

77.1 

159.7 

Optimal 

Design 

0.01 

2.00 

0.01 

9.99 

.0023 

.002 

.0028 

.0014 

.0015 

Table  2.  Program  Data,  Doublet  Input  Response 




NPT 

200 

NA 

5 

NS 

5 

MC 

2 

IA0PT 

1 

INTR 

5 

DELTA 

0.5 

FACTOR 

0 

The  desired  and  final  design  responses  are  given  in  figure  3. 

B.  Optimal  Design  with  Pitch  Command  Step  Input 

In  this  example,  the  control  system  of  figure  2 is  configured  as  follows: 

i)  K and  K.  are  hard-wired  to  zero  ii)  A is  hard-wired  to  1.0.  The  remaining 
z z s ° 

design  parameters  are  chosen  for  optimization.  That  is, 

T 


R = 


k6,kg 


Pitch  Rate  (0) 
Depth  (Z) 


The  command  input  is  taken  to  be  a -25°  degree  pitch  angle  step.  The  desired 
system  responses  are  taken  as  follows 

Pitch  (0)  -25  degree  pitch  angle  (0)  step 

-50  degree/sec.  pitch  rate  (&)  pulse,  one  half 

second  wide  to  allow  leading  edge  of  pitch  angle  step  to  occur. 

Depth  response  (Z)  corresponding  to  the  relationship 

Z = W - U 0 (15) 

o 

Forward  Velocity  u - zero 

Plunge  Velocity  w - zero 

The.  results  of  this  experiment  are  given  in. Table  3;  the  program  settings  are 
listed  in  Table  4. 
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Velocity,  W ‘ I flli  Pitch  Rate 


Comparison  of  Desired  and  Final  Responses,  Feedback  Optimization 


1 


Table  3.  Design  Parameters  and  Errors,  Pitch  Step  Response 


1 

■brfl  Brfl 

1 

RMS  % Difference  (%) 

U W 6 0 Z 

Baseline 

1 

0.0 

0.0 

1.0 

1 

0.0 

0.0 

87.3 

615.5 

Optimum 

Des ign 

5.07 

13.71 

1.0 

100.0 

100.0 

92.1 

5A  2 

106.9 

Table  4.  Program  Data,  Pitch  Step  Response 


Factor 

0 

P 

0 

NS 

5 

MC 

2 

IAOPT  ' 

1 

INTR 

5 

DELTA 

0.5 

KPT 

200 

The  desired  and  final  design  responses  are  shown  in  figure  4.  Note  that  the 
desired  responses  for  0 and  Z are  chosen  so  as  to  be  compatible  v?ith  a -25  degree 
pitch  angle  step.  Observation  of  the  RMSZ  Differences  reveals  that  four  of  the 
trajectories  differ  significantly  from  the  desired  trajectories.  It  should  be 
noted,  however,  that  although  desired  trajectories  were  specified  for  all  five 
states,  only  the  pitch  angle  (0)  was  actually  optimized.  This  is  true  because 
all  entries  of  the  Q matrix  except  Q(4,4),  which  corresponds  to  the  pitch  angle, 
were  hard-wired  to  zero.  The  only  significant  error,  therefore,  is  that  of  the 
pitch  angle  response,  which  is  relatively  small. 

III.  GEOMETRY  AND  FEEDBACK  OPTIMIZATION  (PITCH  MANEUVERS) 

In  this  phase  of  investigation,  a combination  of  feedback  and  geometry 
parameter  optimization  using  the  stern  plane  model  of  figure  2 was  attempted. 

The  experimentation  centers  around  the  development  of  an  effective  method  of 
calculating  the  weighting  matrix  P.  All  case  examples  use  the  pitch  angle 
step  input  and  desired  response  specifications  given  in  the  previous  experiment. 

A . System  Design  witli  P=0 

Let  the  matrix  P be  set  to  zero,  i.e. 


Pfi .il=0  for  all  i 


The  feedback  parameters  and  K.  are  hard-wired  to  zero  so  the  optimization 
procedure  dealt  with  the  following  parameter  set 

R = [K^,K0,?1T  (16) 

The  results  of  the  experiment  are  given  in  Table  5;  the  program  settings  are 
listed  in  Table  6. 

Table  5.  Design  Parameters  and  Errors,  P=0 


Design  Parameters 

Ke  Ko  xs 

RMS  % Difference  (Z) 

U W 0 0 Z 

Baseline 

0.0 

0.0 

1.0 

0.0 

0.0 

87.3 

615.5 

Optimal 

Design 

5327 

3468 

.009 

100 

100 

117.2 

3.5 

106.5 

Table  6.  Program  Settings,  P=0 


FACTOR 

0 

P 

0 

NS 

5 

MC 

2 

1A0PT 

1 

INTR 

5 

DELTA 

0.5 

NPT 

200 

The  desired  and  final  responses  are  given  in  figure  5.  Although  an 
improvement  in  the  pitch  angle  error  has  been  achieved,  the  feedback  parameter 
values  are  quite  large  in  comparison  to  those  of  the  previous  experiment.  The 
program  MOF-NP  has  the  capability  to  penalize  large  departures  from  a set  of 
a priori  parameters.  This  is  achieved  by  calculating  the  weighting  matrix  P of 
equation  (6)  using  the  adaptive  method. 

B . Design  Using  Adaptive-  Method 
Consider  the  matrix  P given  by 


P(i,i) 


FACTOR 
PARER  " 


where 


(17) 


PARER 


NPABC 

I 


k=l 


A A 


2 


(18) 
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The  following  constraints  were  applied  to  PARER  and  FACTOR, 

if  PARER  £ 0.2,  FACTOR  =0.0  (19a) 

if  PARER  ^1.0,  FACTOR  = PARER  (19b) 

if  0.2  < PARER  < 1.0,  FACTOR  = 0.05  (19c) 

The  relationship  between  PARER  and  factor  is  shown  graphically  in  figure  6,  and 
the  relationship  between  PARER  and  P is  given  in  figure  7, 


Figure  6.  PARER  vs  FACTOR 


Figure  7.  PARER  vs  P 


< U°  1.00 

3 0.75 

0 . 50 
£ 0.25 


| h 

0.2  0.5  1.0 

PARER 


It  is  seen,  from  (18),  that  PARER  is  a measure  of  the  normalized  deviations  of 

/*.  >N  A 

c.  from  the  a priori  values,  c . Thus,  without  the  constraints  (19),  as  c 
k /s  k o k 

approached  c , P would  become  large  and  only  small  variations  of  parameters 
from  a priori  values  would  be  allowed.  The  constraints  of  (19),  however, 
frees  the  optimization  process  from  penalties  for  departures  from  a priori 
values  when  the  parameter  estimates  are  close  to  a priori  values,  thus  allowing 
a greater  degree  of  optimization  flexibility. 
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Example  B-l 

Using  the  design  values  of  section  2-B  as  apriori  values, 

= 5.068 

KQ  = 13.713 

X2  = 1.0 
s 

with  K and  hard-wired  to  zero,  feedback  and  geometry  parameter  optimization 
was  performed  for 

A A A A *)  T 

R = lK0*K0’Xs] 

The  results  are  given  in  Table  7;  the  program  setting  are  listed  in  Table  8. 
Table  7.  Design  Parameters  and  Errors,  P adaptive.  Example  B-l 


Design  Parameters 

K0  K0  Xs 

1 

RMS  % Difference  (%) 

U W 6 0 Z 

Baseline 

0.0 

0.0 

1.0 

1 

0.0 

0.0 

87.3 

615.5 

Optimal 

Design 

— 

5.11 

14.41 

1.07 

100 

100 

83.9 

5.11 

106 . 8 

Table  8.  Program  Setting,  Example  B-l 


FACTOR 

0.05 

P 

ADP 

NS 

5 

MC 

2 

IA0PT 

1 

INTR 

5 

DELTA 

0.5 

NPT 

200 

The  desired  and  final  design  responses  are  given  in  figure  8.  A comparison 
between  the  a priori  design  and  the  final  design  is  given  in  Table  9. 


• Table. 

9.  A priori  vs  final  design,  Example  B-l 

Design  Parameters 

Ke  Ke  K 

RMS  % Difference  (%) 

U W 6 0 7. 

A Priori 

13.71 

5.07 

1.0 

100 

100 

92.1 

5.43 

106.9 

Final 

Design 

14.41 

5.11 

1.07 

100 

100 

83.9 

5.11 

106.8 

1 

Example  2 

With  and  hard-wired  to  zero,  feedback  and  geometry  parameter  optimization 
was  performed  for 

* * * '5  T 

R - [KS,K0,Xs2) 

with  a priori 

Ro  = [0,0, 1]T 

The  results  are  given  in  Table  10;  the  program  settings  are  listed  in  Table  11. 

Table  10.  Design  Parameters  and  Errors,  P adaptive.  Example  B-2 


Design  Parameters 

k8  k6  xi 

RMS  % Difference  (%) 

u w e e z 

Baseline 

0 

0 

1.0 

0.0 

1 

0.0  ! 87.3 

1 

615.5 

Optimal 

Design 

0.87 

-0. 11 

1.20 

100 

100  345.9 

16.2 

125.4 

i 

* 

1 


Table  11.  Program  Settings,  Example  B-2 


FACTOR 

0.05 

P 

ADP 

NS 

5 

MC 

2 

IAOPT 

1 

INTR 

5 

DELTA 

0.5 

NPT 

200 

The  desired  and  final  design  responses  are  given  in  figure  9. 

C . Design  using  J0PT2  adjustment 

In  an  attempt  to  improve  the  trajectory  fit  of  the  previous  example,  the 
following  condition  was  imposed, 

if  J0PT2  <_  8,  factor  = 0.0  (20) 

This  condition  completely  frees  the  optimzation  process  from  penalties  for 
departures  from  a priori  values  during  the  first  four  optimization  passes 
(first  eight  iterations).  Using  (20)  in  conjunction  with  the  P adaptive 
method,  (17),  (18)  and  (19),  parameter  optimization  was  performed  for 
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Desired 

Final 


Comparison  of  Desired  and  Final  Responses,  Example  B-2  Adaptive 


r 


with  a priori 

Rq  = [0,0, 1]T 

and  and  K.  hard-wired  to  zero.  The  results  are  given  in  Table  12;  the 
program  settings  are  listed  in  T.'ble  13. 


Table  12.  Design  Parameters  and  Errors,  P adaptive,  J0PT2  adjustment 


1 

Design  Parameters 

Ko  K0  xl 

1 

RMS  % Difference  (%) 

U W 0 6 z 

Baseline 

■ 

0.0 

0.0 

1.0 

1 

0.0 

0.0 

87.3 

r 

615.5 

1 

! 

Optimum 

Design 

1 

9.93 

.115 

7.03 

100 

100 

103.9 

1.19 

105.8  j 

Table  13.  Program  Data,  J0PT2  adjustment 


FACTOR 

0.05 

1 

P 

ADP 

NS 

5 

1 

MC 

2 

IAOPT 

i 

1 

INTR 

5 

DELTA 

0.5 

NPT 

200 

The  desired  and  final  responses  are  given  in  figure  10. 

I 
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IV.  FEEDBACK  AND  GEOMETRY  OPTIMIZATION  (DEPTH  MANEUVERS) 


In  this  phase  of  investigation,  a refinement  of  the  Adaptive  Method 
given  in  section  3.B  is  developed  and  applied  to  parameter  optimization 
for  depth  maneuvers.  The  original  implementation  of  the  Adaptive  Method, 
given  in  (19)  , has  the  upper  bound  for  parameter  variations  controlled  by 

if  PARER  ^ 1.0,  FACTOR  = PARER  (19b) 

where  PARER  is  given  by  (18)  and  P is  given  by  (17).  With  (19b)  in 
effect,  PARER  greater  than  1.0  will  result  in  large  penalties  in  the 
parameter  optimization  scheme.  This  method  is  effective  if  the  priori 
valxies  are  good  estimates  of  the  actual  optimal  system  parameters.  If, 
however,  the  priori  values  are  not  close  to  the  optimal  parameter  values, 

(19b)  may  severely  handicap  Lhe  optimization  procedure. 

Consider,  for  example,  the  case  where  some  parameters  have  priori  values 
equal  to  zero.  In  this  case,  PARER  is  given  by 

-2  (rko  ” rk)2 

PARER  = Er^  + E = PARER]  1-  PARER2  (2]) 

r.  ' 

ko 

where  PARER1  is  determined  by  parameters  with  zero  priori  values  and 
PARER2  is  determined  by  parameters  with  non-zero  priori  values.  If,  at  any 
time  during  the  optimization  procedure,  any  zero  priori  parameter  takes  a 
value  greater  than  1.0  (absolute  value),  (19b)  sets  FACTOR  equal  to  PARER 
and  the  resulting  penalty  is  large.  In  fact,  examination  of  (17)  shows 
that  the  weighting  matrix  P is  no  longer  a function  of  either  FACTOR  or 
PARER  in  this  case,  but  is  set  to  1.0.  This  scheme,  therefore,  will  not 
allow  a zero-priori  parameter  to  exceed  unity  (absolute  value). 

An  improvement  in  the  Adaptive  Method  is  made  when  (19b)  is  changed  to 

If  PARER  2 MAXER,  FACTOR  =1.0  (22) 

where  MAXER  is  a variable,  usually  chosen*  between  1 and  300.  In  (22), 
when  PARER  exceeds  the  designated  upper  bound,  the  penalty  for  departure 
from  priori  values  is  increased  significantly  (factor  typically  is  in- 
creased by  one  to  two  orders  oi  magnitude),  while  the  weighting  P remains 
_ 

MAXER  should  be  ch 'sen  according  to  the  anticipated  variation  of  para- 
meters from  the  sp.  ified  priori  values.  For  example,  a MAXER  of  A00 
would  allow  a zero-priori  parameter  to  assume  values  up  to  20  (absolute 
value) . 
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a function  of  the  parameter  and  priori  values  (PARER).  In  the  case  where 
geometry  parameters  are  specified  with  zero  priori  values,  a MAXER  of  200 
has  been  found  to  yield  reasonable  results  (for  a comparison  of  optimization 
efficiency  of  MAXER  = 1 vs  MAXER  = 200,  see  example  4.B2  and  4.B3). 

A)  Optimal  Design  with  Ramp-Step  Input 

Consider  the  system  given  in  figure  2,  excited  with  the  Depth 
Command  given  in  figure  11  and  specified  as  follows:  (1)  Input  ramp  with 
2.5  feet/second  slope  for  time  ^ 10  seconds,  (2)  Input  constant  at  25  feet 
for  time  >10  seconds. 


30 

4J 

a> 

0) 

20 

4-1 

♦d 

P4 

a 

Q 

10 

t~ 

40 


80  time (sec) 


Figure  11.  Depth  Command  (Zcom) 

In  the  following  examples  the  desired  system  responses  are  taken  a:  follows: 
Depth (Z)  same  as  Depth  Command  (figure  11) 


Depth (Z) 

Pitch(0)  -25/U^  pitch  angle  pulse  for  time  < 10  seconds 

0 elsewhere 

Pitch  Rate(0)  (+)  5.0/uo  pitch  rate  (0)  pulse,  one  half 

second  wide  at  leading  (trailing)  edge  of 

pitch  angle  pulse. 

Forward 

Velocity(u)  zero 

Plunge 

Velocity(w)  zero 

Example  A-l 

In  this  example,  feedback  parameter  optimization  is  performed  for 

R - I'VKo.V'V7 
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with  A hardwired  to  1.0.  The  results  are  given  in  Table  14;  the  program 
s 

settings  are  listed  in  Table  15. 


Table  14.  Feedback  Parameters  and  Errors,  Example  A-l 


RMS  % Difference  (%) 

W 0 0 ' Z 

0.0 

48.38 

74.56 

— 

- 

- 

- 

100 

141.1 

57.67 

4.79 

Feedback  Parameters 


K 


Ka 


Baseline 

0.0 

0.0 

0.0 

0.0 

- 

Priori 

0.0 

0.0 

0.0 

0.0 

- 

1.27 

.029 

1.02 

.079 

- 

Table  15.  Program  Data,  Example  A-l 


FACTOR 

0.01 

1 

P 

ADP1 

NS 

5 

I 

MAXER 

200 

1A0PT 

1 

I 

MC  ’ 

2 

DELTA 

0.5 

I 

NPT 

200 

MAXR 

16 

1 

UNKR 

m 

The  desired  and  final  responses  are  given  in  figure  12. 

Example  A-2 

In  this  example,  feedback  and  geometry  optimization  is  performed  for 

A A A A A ^ T 

R = fKz,K0,KZ’K0,As^ 


20  - 


The  results  are  given  in  Table  16;  the  program  settings  are  listed  in 
Table  17. 


Table  16.  Design  Parameters  and  Errors,  Example  A-2 


Design  Parameters 

h h H H xl 

KMS  % Difference  (%) 

W 6 0 Z 

Baseline 

0.0 

0.0 

0.0 

0.0 

1.0 

0.0 

48. 38 

74.56 

— 

Priori 

0.0 

0.0 

0.0 

0.0 

1.0 

- 

- 

- 

- 

Optimal 

Design 

.093 

.027 

.214 

.131 

1.76 

100 

137.9 

57.8 

3.94 

Table  17.  Program  Data,  Example  A-2 


— 

FACTOR 

0.01 

P 

ADP1 

NS 

5 

MAXER 

200 

IAOPT 

1 

MC 

2 

DELTA 

0.5 

NPT 

200 

MAXR 

16 

UNKR 

5 

The  desired  and  final  responses  are  given  in  figure  13.  Notice  that  tin 
optimal  values  obtained  for  K^,  K*  and  look  reasonable,  but  seems 

inappropriately  small.  This  is,  however,  a depth  manuever,  and  as  such  the 
parameters  obtained  are  acceptable,  although  it  is  possible  that  the  pitch 
response  with  this  set  of  parameters  is  poor. 


t 


Example  A-3 


In  this  example,  feedback  and  geometry  optimization  is  performed  for 

A A A A A A n 'T* 

R - [K2,K(|,KZ,K0,X2] 

where  the  priori  values  for  the  parameters  are  based  on  the  results  obtained 
in  section  3C(pg.  16).  The  results  are  given  in  Table  18;  the  program 
settings  are  listed  in  Table  19. 


Table  18.  Design  Parameters  and  Errors,  Example  A-3 


Design  Parameters 

Ke  Kz  k6  Ki  xl 

RMS  % Difference  (%) 

W 0 0 Z 

Baseline 

9.93 

0.0 

0.115 

0.5 

7.03 

0.0 

48.38 

74.56 

— 

Priori 

9.93 

0.0 

0.115 

0.5 

7.03 

- 

- 

- 

- 

Optimal 

Design 

9.07 

.213 

.115 

.442 

6.57 

100.0 

93.2 

61.8 

2.83 

Table  19.  Program  Data,  Example  A-3 


FACTOR 

0.01 

ADP1 

NS 

5 

MAXER 

200 

IAOPT 

1 

MC 

2 

DELTA 

0.5 

NPT 

200 

MAXR 

16 

UNK.R 

5 

The  desired  and  final  responses  are  given  in  figure  14.  Notice  that  the  values 

2 

obtained  here  for  Kq  and  are  quite  different  from  those  obtained  in  the 
previous  example.  These  parameters  are  comparable  to  those  obtained  in  the 
pitch  optimization,  and  therefore  may  be  acceptable  for  both  depth  and  pitch 
maneuvers . 
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B)  Optimal  Design  with  Triplet  Input 


In  this  phase  of  investigation,  the  stern  plane  feedback  system 
is  excited  with  the  depth  command  (Zc  ) input  given  in  figure  16.  The 
desired  depth  response  is  the  same  as  the  depth  command. 


Figure  16.  Depth  Command  Input  (Z  ) 

com 


Example  B-l 

In  this  example,  feedback  and  geometry  optimization  is  performed  for 


K = [K*,Ka,K,,K0,A 


2 ,T 


'z’*6‘  z*  oi  s- 

with  the  feedback  priori  values  set  to  zero.  The  results  arc  given  in  Table  20; 
the  program  settings  are  listed  in  Table  21. 
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Table  20.  Design  Parameters  and  Errors,  Example  B-l 


Design  Parameters 

Ke  Kz  H Ki  K2 

Baseline 

0.0 

0.0 

0.0 

0.0 

H 

Priori 

0.0 

0.0 

0.0 

0.0 

^oj 

Optimal 

Design 

.371 

.074 

.296 

.106 

RMS  % Difference  (%) 

W 0 0 Z 

0.0 

96.7 

186.6 

— 

[■l 

- 

- 

- 

100 

152.9 

59.5 

23.63 

Table  21.  Program  Data,  Example  B-l 


FACTOR 

0.01 

P 

ADP1 

NS 

5 

MAXER 

1 

IAOl’T 

1 

MC 

2 

DELTA 

0.5 

NPT 

200 

MAXR 

16 

UNKR 

5 

The  desired  and  final  response  are  given  in  figure  17.  Notice  that,  in  this 
example,  MAXF.R  equals  1 and  all  of  the  final  feedback  parameter  values  are  less 
than  1.0. 

Example  B- 2 

This  example  is  identical  to  the  previous  example,  with  the  exception 
that  MAXER  is  set  at  200  instead  of  1.  The  results,  given  in  Table  22,  show 
a definite  improvement  in  the  depth  response  error.  Notice  that  some  of  the 
feedback  parameters  (Kg,K^)  have  absolute  values  greater  than  1.0. 
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Desired 

Final 


Comparison  of  Desired  and  Final  Responses,  Example  B-l 


Table  22.  Design  Parameters  and  Errors,  Example  B-2 


Design  Parameters 

K0  Kz  H Kj  X* 

RMS  % Difference  (%) 

W 9 0 Z 

Baseline 

0.0 

0.0 

0.0 

0.0 

1.0 

0.0 

96.7 

186.6 

— 

Priori 

0.0 

0.0 

0.0 

0.0 

1.0 

- 

- 

- 

■EZSfSI 

H 

.085 

-4.12 

■9 

6.51 

100 

131.1 

43.5 

11.06 

Table  23.  Program  Data,  Example  B-2 


FACTOR 

0.01 

1 

P 

ADP1 

NS 

5 

I 

MAXER 

200 

IAOPT 

■1 

I 

MC 

2 

DELTA 

0.5 

I 

NPT 

200 

MAXK 

16 

1 

UNKR 

5 

The  desired  and  final  responses  are  given  in  figure  18. 

Example  B-3 

In  this  example,  feedback  and  geometry  optimization  is  performed  for 
R - [K^,K^,Kz,Kq,A^]T 

with  priori  values  for  the  parameters  based  on  the  results  obtained  in 
section  3C  (pg.  16).  The  results  are  given  in  Table  24;  the  program 
settings  are  listed  in  Table  25. 
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Table  24.  Design  Parameters  and  Errors,  Example  B-3 


Design  Parameters 

Ke  S K6  Kz 

Baseline 

9.93 

0.0 

.115 

0.5 

2.03 

Priori 

9.93 

0.0 

.115 

0.5 

7.03 

Optimal 

Design 

5.48 

1.106 

.113 

.224 

9.88 

RMS  % Difference  (%) 

W 0 0 Z 

0.0 

96.7 

186.6 

B 

■ 

- 

- 

- 

100 

105.4 

46.3 

6.82 

Table  25.  Program  Data,  Example  B-3 


EACTOR 

0.05 

P 

ADP1 

NS 

5 

MAXER 

200 

IAOPT 

1 

MC 

2 

DELTA 

0.5 

NPT 

200 

MAXR 

16 

UNKR 

5 

The  desired  and  final  responses  are  given  in  figure  19. 

I 
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V. 


PARAMETER  OPTIMIZATION  UTILIZING  EXPONENTIAL  TRANS FORMAT  ON 

In  this  phase  of  investigation,  an  exponential  transformation  is  implemented 

to  guarantee  the  acquisition  of  non-negative  design  parameters.  Consider  the 

transformation 

~ rj  r.-  • r.j 
p . = e -L  = c 1 x 


9Pl  _ 


9ri 


= ri  e 


(23) 


where  p represents  the  actual  physical  parameters  and  r represents  the  transformed 

variables.  This  transformation  provides  an  isomorphic  mapping  from  the  real 

numbers  (r)  to  the  non-negative  real  numbers  (p).  Thus,  if  optimization  is  now 
T 

performed  on  R , the  optimal  design  parameters  must  necessarily  be  non-negative. 

It  should  be  noted  that  each  element  in  t lie  Ii?  matrix  given  in  (10)  can  be 
calculated  as  follows 


9c  _ 9c  * ?[j>  _ 9c 
9?  9p  9ir  = 9p 


(24) 


A difficulty  encountered  in  implementing  this  transformation  is  the  selection 
of  limits  for  allowable  values  for  r and  p.  The  necessity  of  limits  arises  from 
the  following  observations:  (1)  as  p approaches  zero,  r approaches  negative 
infinity,  (2)  if  p takes  on  values  close  to  1.0,  r becomes  very  small, 

(3)  moderate  values  of  r greater  than  1.0  will  result  in  large  values  of  p. 

The  limits  chosen  to  control  the  parameters  are 


if  r 

< -18.4  , 

r = 

-18.4 

(25a) 

♦ 

if  1 

n 

A 

O 

1 

— i 

o 

r = 

0.0 

(25b) 

if  r 

>6.0  , 

r = 

6.0 

(25c) 

Equations 

(25a)  and 

(25c) 

define  the  range  of  allowable  r values. 

This 

I 

corresponds  to 

p 

'max 

- 403 

'min 


10 


Equation  (25b)  restraints p from  assuring  values  very  close  to  1.0,  while 
still  allowing  it  to  be  exactly  1.0. 
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Examp le  1 


r 


In  this  example,  feedback  and  geometry  optimization  is  performed  for 

T 


where  the  depth  command  (Z  ) and  desired  depth  response  is  given  in 

com 

figure  16.  The  results  are  given  in  Table  26;  the  program  settings  are  listed 
in  Table  27. 


Table  26.  Design  Parameters  and  Errors,  Example  1 


Design  Parameters 

Ko  Kz  Ko  xs 

RMS  % Difference  (%) 

W 0 0 Z 

Baseline 

0.01 

0.01 

0.01 

0.01 

1.0 

0.0 

20.99 

99.43 

98.83 

Priori 

0.01 

0.01 

0.01 

0.01 

1.0 

- 

- 

- 

- 

Optimal 

Design 

9.38  . 

1.31 

0.01 

.0068 

7.43 

100.0 

166.0 

49.14 

6.62 

Table  27.  Program  Data,  Example  1 


FACTOR 

0.01 

P 

ADP1 

NS 

5 

MAXER 

10 

IAOPT 

1 

MC 

2 

DELTA 

0.5 

NPT 

200 

MAXR 

16 

UNKR 

5 

Notice  the  selection  of  Baseline  and  Priori  parameters  were  chosen  to 
avoid  the  problem  areas  controlled  by  the  limits  in  (25).  The  one  exception 
is  the  selection  of  1.0  for  \ * which  yields  and  r value  of  exactly  zero.  The 
desired  and  final  responses  are  given  in  figure  20. 


Desired 

Final 


Compairson  of  Desired  and  Final  Responses,  Example 


Example  2 


This  example  is  identical  to  the  previous  example  with  the  exception  that 
2 

the  prior  guess  for  is  chosen  to  avoid  1.0.  The  results  are  given  in  Table 
28,  the  program  settings  are  listed  in  Table  29. 

Table  28.  Design  Parameters  and  Errors,  Example  2 


Design  Parameters 

Ke  Kz  Kg  KJ  X2 

RMS  2 Difference  (2) 

W 0 0 Z 

Baseline 

0.01 

0.01 

0.01 

0.01 

1.0 

0.0 

20.99 

99.43 

98.83 

Priori 

0.01 

0.01 

0.01 

0.01 

0.1 

- 

- 

- 

- 

Optimal 

Design 

16.36 

2.82 

.0091 

.0072 

17.63 

100 

158.9 

40.92 

5.02 

Table  28.  Program 


FACTOR 

0.01 

NS 

5 

IAOPT 

1 

DELTA 

0.5 

MAXR 

16 

Data,  Example  2 


P 

ADP1 

MAXER 

10 

MC 

2 

NPT 

200 

UNKR 

5 

Note  that  while  the  values  of  K*  and  K’  are  practically  the  same  as  in  the 

2 

previous  example,  the  values  of  Kg,  K^,  and  Ag  increased  significantly,  resulting 
in  a slightly  better  depth  response.  The  desired  and  final  responses  are  given 
in  figure  21. 
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Desired 


Figure  21.  Comparison  of  Desired  and  Final  Responses,  Example 
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APPKNDIX  A 


Description  of  Design  Parameters 


Geometry  Parameters 


Coning  tower  height 


Bow-plane 


Stern-plane 


Rudder 


Stern  Plane 


R.  , = K 

z 


*13  “ K1 


APPENDIX  B 

Iteration  Data  for  MOF-NP 


KOPT 

NTER 

JOPT1 

JOPT2 

INTR 

IPROPT 

IAOPT 

1 

2 

0 

1 

1 

1 

1 

0 

2 

0 

2 

1 

1 

1 

0 

2 

0 

4 

1 

1 

1 

0 

2 

0 

8 

1 

1 

1 

0 

2 

0 

16 

1 

1 

1 

0 

2 

0 

32 

1 

1 

1 

0 

2 

0 

64 

1 

1 

1 

0 

2 

0 

128 

1 

1 

1 

0 

2 

0 

500 

1 

1 

1 

0 

2 

1 

0 

1 

1 

1 

0 

2 

1 

0 

2 

1 

1 

0 

2 

1 

0 

4 

1 

1 

0 

2 

1 

0 

5 

1 

1 

r 


APPENDIX  C 

The  linearized  state  equations  of  a vehicle  are  of  the  form 
dx 

A^p-  = Bx  + Cu 
at 


where 


u = D 


com 


+ Ex 


Equation  (C.l)  can  therefore  be  rewritten  as 

A~ ■ = (B+CE)  x + CD  „C°m 
at  Z 

com 

The  matricies  of  (C.3)  are  as  follows: 


0 

0 

0 

0 


w 


0 

m-Z 

-Mw 

0 

0 


0 

-Z- 


1 -M- 

y q 


o 

o 


z +mU 

q o 


1 

0 


- A3  - 


0 

0 

0 

1 

0 


0 

0 

0 

0 

1 

0 

0 

0 

0 

0 


(C.l) 


(C.  2) 


(c.3) 


(C.4) 


(C . 5) 


0 

0 


0 

0 


1)  = 


-k; 


-Ke 


K.S 


E 


0 

0 


-K0b-  K^Uo 


Ko  + KzSl,» 


Next,  redefine  matricies  >’  and  C as  follows; 


bnew  = B+CK 


CNEW  = CD 


Thus  (C.3)  becomes 


Adt  BNEW  X + ( NEW 


com 


com 


44  - 


K., 


-K, 


(C.6) 


(C . 7) 


(C.8) 


(C.9) 


(C.10) 


1 


Therefore  we  have 


r 


APPENDIX  I) 

The  longitudinal  dynamics  of  the  USl'-KPV  vehicle  are  governed  by  the 
vector  state  equation 


r ■ 


where 


m 

I 


Z. 

w 


M. 

q 


vx-  z . 
w 

-z. 

q 

0.0 

0.0 

W 

i -M. 
w 

l -M. 

y q 

0.0 

0.0 

_d_ 

6 

-0.0 

0.0 

1.0 

0.0 

dt 

e 

0.0 

0.0 

0.0 

1.0 

z 

z 

w 

Z +mU 

q o 

CM 

CD 

0.0 

w 

X 

Z6 

s 

'6h' 

M 

w 

M 

q 

Me 

0.0 

6 

+ 

MS, 

b 

M6 

s 

5 

t sj 

0.0 

1.0 

0.0 

0.0 

0 

0.0 

0.0 

1.0 

0.0 

-u 

o 

0.0 

z 

_ 

! 

0.0 

0.0 

= 280.39  slugs 

= 14267.0  slug  ft2 

= 14.616  ft/sec. 

= -270.87369 

= M.  = -49.718618 
w 

= -12818.43975 

= -[ (5.693230U  ) + X2  (.037748  U ) + A2  (4.76436  U ) 

O C Ob  o 

+ X2  (6.01750  U ) + x2  (.030661  U )] 

S O K o 

= (192.15600  U ) + A2  (23.095907  U ) + A2  (-85.5179  U ) 

o B o s o 

= — l (56.4  697  U ) + A (-22.26927  U ) + \2  (11 1 .(>.'7245  r ) ] 

o B o s o 

(-931.456809  U ) + A2  (-108.451145  U ) + A2  (-1596.298837  U ) 
o B os  o 

- 0.0 

= -722.2 

- -A2  (5.591056)  I’2 

b o 
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APPENDIX  D 
- Continued  - 


Z6  = -X2  (2.503594)  U2 

s s o 

M6b  = A2  (27.228202)  U2 

M6  = A2  (-37.203137)  U2 

s s o 


I 


r 


0011  T£ll=l  ,NTER 
HRITE(u,403)  ITER 

FORMA  f(/lX,20( SUITER  ) , /SOX,  1 I TFRAT  I OM  NO. MS, 
l/3JX,lo(lil-),/50X,13(  111-)) 

I F ( ITRAU.EQ.  J)  CO  TO  SO'J 
CO  iJJ  1=1, IUNKR 
J=INDXR(  I ) 

IF(R(J).GE.u.O)  R(J)=b.O 

I F(R( J) . L£ . -13.4)  R(J)=-18.4  THIS  PAC 

IF(RC(J) .GE.b.O)  RC(  J ) =b . 0 , ^rvv 

I F( RC( J ) . LE . - IS . 4 ) RC ( J ) =-  13 . 4 FROM 

IF(AD3(iUJ)).LT.l.E-lJ)  R(j)=0.0 

COOT  I RUE 

CONTINUE 

CALL  SAVER 

I LOG2=0 

I P0S=1 

CALL  SELECT 

CALL  ERROR! YY, SUMER,  I Oil) 

CALCULATE  THE  .iEIOilTIIIG  MATRIX  Q USED  IN  J1 

001 Jj I =1,0  TM 
TEM“QQ( I ) 

I FCTtM. EQ. 0 . )TEM=1 . . 

Q( I ) =1 . / ( FLOAT (NTH) *TEM) 

ac i ) =o . o 

Q( 2) =0 . 0 
0(3) =0.0 
Q( 4) =0 . 0 
Q( 3}=0 . 0 

HRITE(b,  130G)Q(  I ) 

CONTINUE 

FORMAT! 10X, 'Q  EQUALS  ’,G14.6) 


»IsPAOEisBtsTW«.tnra»cnc*«u 

FROM  COr'X  FURMISH®  10  DOQ  " 


CALCULATE  THE  HEIGHT  I NG  MATRIX  Q2  USED  IN  J2 

PARER=0 . J 

ZAUER=J.O 

UARER=0 . 0 

DOS  03 I PMC-1, IUNKR 

J=INDXR( I PMC) 

I EM* ( PR  I ORE ( I PMC ) -R( d ) *RC ( 0 ) ) » * 2 
TEM1=(PRI0KE( IPMC))**2 
I F ( TEM1 . NE . 0 . 0 )GO  TO  DOS 
TEi  11  = 1.0 

ZARER’ZAREK+TEM/TEMl 

CONTINUE 

Q2 ( I PMC ) =1 . 0/TEM1 
PARE  R = PARE  KfTEM/ TEill 
HARE R= PARE R-2ARER 
IF(IADPT.EQ.O)  GO  TO  Db2 

ADAPTIVE  METHOD 

RARER=PARER 

FACT0R=TEM3 

PM AX  IS  PROGRAM  VARIABLE  FOR  MAXER 
PIIAX  = 200.0 

I F ( I TRAN. EQ. 1)  PMAX  = 10.0 
I F (RARER. LL. 0. 2) FACTOR’D. 0 
I F ( RARER . GE . PMAX ) FACTOR’l . 0 
I F ( TEMD . LE . J . 00001) FACT0R=TEM8 


THIS  PAGE  IS  BEST  QUALITY  PRACTICABLE 
FROM  COPY  FUWOLUED  TO  LOG - 


u390  9o2  CONTINUE 

o90J  I F ( I PUNT . AGE  . 1 ) NR  I TE( 6, 930) PARER 

o J1U  I F(  I PRNT.GE.  1)17RI  TE ( 6, 966 ) 2ARER 

U J2  J I F ( I PRUT  .GE.  1 )NRI  TE(6,967  ) OARER 

0930  930  FOUMAT ( 1 7X,  1 PARER  = ',0.14.0) 

0940  900  FORMAT(17X,  ' PARER1= ' ,014 . 6 ) 

O3:>0  907  FOUMAT (17X, 'PANE I! 2=' ,014.0) 

OOuO  IF(PARER.EQ.O.O)  PAUER'F  LOAT  ( I UNKIi) 

u370  DOSOulPMC-l, IUNKR 

0930  IF(IAOPT.LE.l)  00  TO  963 

09 3 J C 

7000  C JOPT2  ADJUSTMENT 

7010  C 

7020  I F ( JOPT1 . EQ. 1)00  TO  960 

iojO  IF(JOPT2.LE.3)FACTOR*0.0 

7040  900  CONTINUE 

70^0  3o3  CONTINUE 

70u0  C 

70  70  02(1  PMC ) =Q2( I PMC) * FACTOR/ PAREU 

7030  i/RITE(o,1901 ) I PMC,Q2(  I PMC) 

709 J SJu  CONTINUE 

7100  1301  FOUMAT ( 20X, 12,'  Q2  EOUALS  ’,014.6) 
7110  U04I-1, IMP) 

7120  D04 J=1 , NTP1 

7130  4 G( I , J ) =0 . 

7140  DOS  I =1,NT 

7130  X(l)=0. 

7100  'j  XD(  I )=0. 

7170  DOo I =1, NA 


THIS  PAGE  IS  BEST  QUALITY  PRACTICABLE 
mM  COPY  FURNISHED  TO  DOC  ^ 


1 J34J 
10350  C 
1 J 3 o 0 
1J370 
lJjdd 
1 J5JJ 
10400 
I J **  1 J 
1J42J 
1 JUDO 
1 9 4 4 J 
l J*  > J 
I J4oJ 
l-K/J 
lj^tg  J 
1J-.JJ  c 
lJjdd  c 
I J j I J c 
lja2d  C 
1 J j j J c 
1 J j4  J C 
iJ'jjJ  C 
I J j u J 
1 J > / J 
1 J 3d  J 
ljjjj 

1 Ju  JO 

ijuld 
1 Ju2  J 

l Ju  Jd 

1 J U ‘f  J J 4 'J 
1 Ju  3 J 
1 Joud 
I Ju70 
1 J od  0 J4u 
1 Ju  JO  94  7 
1 J/Jd 
1J/1J 
1 J 720 
1 J 7 3 J C 
i J /**  J C 
1j/jJ  c 
U/ud  C 
i j/7j 
I J7dd 
1 J7J  J 

ij:jj  c 

Udld  C 
ldd  2d 
l'Jgid 
1 Jd-id 
lJdjJ  C 
1 J u u d 
iJd/d 
1 J u d d C 
1 Jd  Jd 
1 J J Jd 
1 J lid 
1 Jj2d 
IjJjd  1 
i j j 4 J C 
lJjjd  C 
I J ju  J 
1 JJ/d 
1 Jldd 
1 JJJd 


SUBROUTINE  PDF.RCR 
IMPLICIT  REAL*3( A-H,0-2) 

C0MM0N/W0RK4/  IUIIKR,  I POS,  I TUAN 
C0MM0N/U0KK5/ INDXIUIO) 

COMMON/NORK3/R(lb),RC(  1G)  ,ll(  2 0, 20  ) , RS(  lb)  ,IIT<  20, 20)  ,GN(  20, 20)  , 

IP( iu),PC(lb),PS(lb) 

COMMON  /MATU IX/A(ld,IO),AI ( 10, 10) , B( 10, 10) ,C( 10,15), 

1AS( Id, Id), 05(10,10), 01 (10, 10), 02(10, 10), CS(10, 15), 01(10,15), 
2C2(ld,la),CA(ld,ld),C0(10, Id) ,CC( 1 J, 1 5 > , X I NT( 10 ) , SB(  10),S03(10), 
jdlNP(ld),XINTo(id),S0S(ld),SI)0S(ld),0l,JP5(ld),Q(20),QQ(20),Ql(20) 
COMMON  / INTUS/MAX,NA,MC,NPUTS,  I PAI!M(  100, 5 ) ,NPAB,  NPABC,  NTER,  I NIK, 
idOPTL, JOPT2,  I AOPT,MAXG, I S( 10) , »S0( Id), I 5M(10), NS, I5DIK  Id) ,MSD, 
2INTS(1G),NI  , I3U(  Id),  MSI),  I SDU(  10  ) , NSOd,  I NPO  ( lu  ) , N I NB,  NPT,MAXNPT, 
3NTP, N TP1, NT, NTM, N 31, N 32 , NS3, NS4,  I PUNT,  I L001,  IL002,  IADPT 
0 I MENS  I ON  I NU(  Id ) ,113  ( 20, 20 ) 


Til  1 5 SUBROUTINE  CALCULATES  THE  PARTIAL  DERIVATIVE  OP  C 
(THE  VECTOR  OP  UNKNOWN  PARAMETERS)  WITH  ULSPECT  TO  R 


DO  045  11=1, IUNKK 
I * I NDXIU  I I ) 

RRC=K( I )‘UC(  I ) 

I F( RRC.GE . o . J ) RRC=b.d 
P( I ) =EXP(RRC) 

PC ( I ) = 1 . d 

IF(ITRAN.EQ.O)  P(I)=R(I) 

I Ft  I TRAN.  El).  0)  PC  ( I )=RC(  I ) 

CONTINUE 

NRI TE( u, 04u ) 

CALL  PRVECt  P, lb ) 

WRI TE ( b, 947 ) 

FORMAT ( /10X, 'P  VECTOR',/) 
FORMAT(/ldX, ' R VECTOR',/) 
CALL  PKVEC ( R, lb ) 

UC0M=14 . u lo 
UCOM2=UCOM«*2 


FORM  X DELTA13  ETC. 

XDI3  = -P(2)*PC(  2)*.  5591db*UCOM2 
ZOIi--P(E)*  PC  (2)*5. 591db  *UCOM2 
MDB  = P( 2)*PC(2)*27.223202  *UCOM2 


XD3=-P(3)*PC(5)*.250359*UCOM2 
ZD5--P( 3)*PC(3)*2.5d3594  *UCOM2 
MUS-P( j)*PC(j)*(-3/.20IlS7) *UCOM2 

TT=-P(13)*PC(13)*P(ll)*PC(ll) 

S3=P( 1 4 ) *PC ( U)*H(12)*PC(12)*P( J)*PC( J)*UCOM 

DO  11*1, 2d. 

DO  10  * 1 , NPAUC 
IK  I , d ) *0 . 0 
H3( I , d ) »d . 0 
CONTINUE 


I ND ( 1 ) * 1 
I ND( 2)*2 
I ND( 3 ) *3 
I ND( 4) *4 
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i 


A 


20000 
2 J JiJ 
2 d d 2 J 
20350 
2 d d 4 d 
2 JiJ  Jd 
2 0 0 b 0 

2dd/d  C 
2 ddd d C 
2dd  Jd 
2 d 1 J J 
2dlld 
2 d 1 2 d C 
2d  1 30  C 
2di*id  C 
2d  lad 
2 d 1 u d 
2dl/d 
2 d i o d C 
2diJd  C 
2J2J0 
202id 
20220 
2023d  C 
20240  C 
2 0 230 
202od 
20  2 70 

2 d 2o  0 
2d2  JO 
20300 
20310  C 
2J320  C 
2 0 330 
20340 
2 J j j0 
203od 
2 J 3 /0 
20300 
2b3  Jd 
2 0 4 d 0 
2 34  Id 
20*20 
2043  J 
20440 

2 0 4 3 0 

2 d ‘V  o o C 
20470  C 
2 040  0 
204  JO 
20sdd 
20310 
20  j20 

2 0330 
20340  C 
20330  C 
2 0 3 o 0 
2 03  7 0 Jlj 
2 J 3 o 0 
2 J 3 J d 320 
2duJ0  JIJ 
2doi0 
2do2  J 
2 Ju3d 
2do40 
2do3d  023 


11Jl)(,,)=l)  THIS  PAGE  IS  BEST  QUALITY  FRACII CABLE 

H(1,1)=-PC(9)*XI)S  FROM  OOfi  FURNISHED  TOLDC  

11(1,3)=  PC(9)*UCOM*XDS 
tt(  1,J)=-PC(9>*2DS 
(1(1,/)=  PC(‘J)*UCOM*ZDS 
H(1,9)=-PC(9)*MU$ 

11(1,11)=  PC  ( 0 ) * UC0M*MQS 


11(2, 2)=  PC  ( 10  ) *XDS 
11(2,0)=  PC  ( 10  ) *ZDS 
11(2,10)=  PC  ( 10  ) *l-iDS 


11(3,4)=- PC (11)  *XDS*P(  13)  *PC(  13 ) 
ll(3,0)=-PC(  11)  »ZOS*P(  13)  *PC(  13) 
11(3, 12)=- PC  (11)  *MDS*P(  13)*  PC  (13) 


H(3, 14)=  PC( 1 1 ) *XDS 
11(3,  lo)=  PC ( 1 1 ) *2DS 
11(3,10)=  PC  ( 1 1 ) *MDS 


11(4,3)=  PC  ( 12)  *XDS*P(  14 ) *PC  ( 14) 
11(4,7)=  PC(12)*Z03*P(14)*PC(14) 
11(4,11)=  PC(  1 2 ) *KD  S * P ( 1 4 ) * PC ( 1 4 ) 
H(4,13)=-PC(12)*X1)S 
H(4,15)=-PC(12)*ZDS 
H(4, 1 / ) = — PC (12) *MDS 


IF(Ai5S(P(3)).LE. 0.00001) GO  TO  913 
H(3,1)=-P(U)*PC(9)*XDS/P(3) 

IK  J,2)=-0.01*PC(3)*(lll.l>2724S*UCOM)+P(10)*PC(10)*XDS/P(3) 
H( 3, 3) =SS*XDS/P( 3 ) 

IU 3,4) =TT*XOo/ PC  3 ) 

ll(3,5)=-1.0*PC(3)*(6.01750*UCOH)-P(9)*PC(9)*ZDS/P(3) 
ll(5,b)=- PC (3)*(111. 027243* UCOM)+P(10)» PC (10)*ZDS/P(3) 

H( 3, 7 )=SS*ZDS/ P( 5) 

H(3,3)=TT*20S/P(3) 

11(3,9)=  PC  (3)*(  -35 . S179*UCOM)-P(  9 ) *PC(  9 ) *MDS/P(  3 ) 

11(3,  10 )=PC(5)*( -1590. 293837 *UCOM)+P( 10)  *PC(  10)  *MDS/P(  3 ) 

H( 3, 11 ) =5S*MDS/ PC 3 ) 

11(3,12)  =TT*WDS/P(  3 ) 


H(5,15)=-P(12)*PC(12)*XDS/P(3) 
11(3,1 4 )=P(11)*PC(11)*XD3/P(3) 
H(3,13)=-P(12)*PC(12)*ZOS/P(3) 
H(5,lb)=P(li)*PC(ll)*ZDS/P(i) 
IK  3, 17 ) =-P( 12 ) *PC( 12) *MDS/ P( 3 ) 
1I(3,1o)=P(11)*PC(11)*MDS/P(3) 


GO  TO  JIJ 
Cd.'ITI  HUE 


>JK  ITE(u,920) 

FOUKAT  ( / 10X,  '**«*«*** ****ERROR 
CONTINUE 

DO  925  1=1, IUNKR 

00  925  J = 1,11PA0C 

1 I = I NO ( I ) 

H3( I ,0 ) =H( I l,J) 

CONTI  HUE 


IN  PDERCR-P(3)=0*********** 


/) 
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1 


CTISPASSISBm3^«™^* 

jHO*  EUEOUBHED  10  


20uu0  00  92b  1=1,5  . „„.,, 

20o/0  DO  92b  J = 1 , NPAIIC  THIS  PAGE  IS  BEST  QUAL1 

20ud0  ll(l,J)=d.U  pt-rtu  f-jOPY  FUKiilSiiED  TO 

2doJ0  92b  COOT  I NUC 

20/00  DO  02/  1=1,1  UiJKI! 

20/10  DO  027  J=1,NPABC 

20/20  il(  I , J)=,I3(  I ,J) 

20/00  02/  CONTINUE 
20/**0  Will  TE  ( o,  OH  1 ) 

2 0 / j 0 041  FORMAT (// 20X,  1 THE  II  MATRIX  AFTER  ili  ADJUSTMENT',/) 

2 07ol)  DO  04  2 l=l,IUNf(R 

20/70  WRI TE ( b,  043)  GK  I , J ) , J = l,  IJPABC  ) 

20/ o0  042  COilTlilUE 

20/00  040  FORMAT ( IX, 3014.0) 

2 Oo 0 0 DO  044  1 =1,  1 UiJKR 

20^10  11=1 NDXR( I ) 

2o0 2 0 PRC=P( I I )*PC< I I )*RC( I I ) 

2 Ou  0 J IF(ITRAN.EQ.O)  PRC  = 1.0 

20040  DO  044  J = l,ilPAUC 

20u  j 0 II C I , J ) =U(  I , J ) * PRC 

iOuuJ  044  COilTlilUE 
20o  / 0 .IRlTE(b,948) 

20u o 0 DO  040  1=1,  IUNKR 

200  00  UR  I TE  ( u,  OoO ) (ll(  I , J ),  J = 1,  NPA3C ) 

20  000  9*9  CONTINUE  . 

2O010  040  FORMAT(/20X, 'THE  115  MATRIX  AFTER  TRANSFORMATION'  /) 

2 0 0 2 J JjJ  FORMAT ( IX, 8G14 . b) 

2000 J RETURN 

2d  040  END 


2JJjJ 
2 0 J o 0 C 
2JJ/0 
2JJ0  0 
2JJ  JO 
2 1 j J J 
2 1 j 1 J 

21 J2  J 
21 J j J 
21040 
21050 
21  JoO 
ZiJ/j 
21  030 
21JJJ 
211J0  C 
21110  C 
21120  C 
2 1 1 j 0 C 

2iuo  c 
211)  0 C 
2 1 ioO 
21170 
21100 
211  JO 
2 120  J 
21210 
21220 
21210 
21200 
21210 
2 1 2 w J ■* 
212/0 
2h.J  J 
2T2  JO 
2 1 j 0 J 
2111 J 
21)20 
21)10 
21)tO  10 

2 J >0  11 
2 1 3oO 
21170 
2 1 ) ) 0 12 
21)0  11 
t »*tOO 
21 .10 
21-i20 

21410 

2 1 4 4 0 l 
210)0 
2 l4oO 
214/j 
21000  2 
210  JO 
21)00 
21)10 
21)20 
2 1)  iO 

2 1 ) 4 0 
21))0 
2 l)uO 
2 1 ) / 0 
21)00 
2 1 1 JO 
2 loJO 
2 lo  1 J 
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SUBROUTINE  DCS  I ON 
IMPLICIT  REAL*0  ( A-ll,  0-  2 ) 

COMMON/ NOIUd/R(  ill),  UC  (10), 11(20,20)  / KSllO/HT  (20, 20),  GN  (20,20), 
lP(lo),PCdb),PS(16) 

COMMON  /MATK I X/A(  10,  10  ) , A I ( 10,  10 ) , B(  10,16 ) , C ( 10, 15  ) , 

1 A 0(10, 10), 00(10, 10), 01(10,10),  I!  2(10, 10), CS(10, 15), Cl(10, 15), 

2C2( 10, 15),CA(10,1U),CB( 10,10),CC(10, 15 ) , X I NT< 10 ) , SB( 10 ) , SDB( 10 ) , 
jOINP(iJ),XIMTS(10),SOS(  10),SD0S( 10 ), 0 1 UPS ( 10 ), Q( 20) , QO( 20 ), Q1 ( 20) 
COMMON  / I NTGS/MAX, MA,HC,NPUTS, I PARI-U  1 00, 5 ) , IIPAC,  NPAUC,  IITER, I MTU, 
IJOI’1'1,  JUPT2,  IAUPT,MAXG,  10(10),  130(10),  I5M(10),NS,  I SUM( 10 ) , N5D, 

2 I .JTS(  10 ) , ;j  I , ISU(10),NS0,  I SOL!  ( 10  ) , NSDO,  INPU(10),.Mi;i0,NPT,MAXNPT, 
3NTP,NTP1,NT,NTM,N31,N32,NS  j,;JS4,  I PRNT,  I LOG1,  I LOG2,  IADPT 
COMMON/ KF0U/E(5,5),CE( 5, 3), 0(5,5) 

COMMON/ NORK4/  I UNKR,  IPOS,  I TRAN 
C0MM0N/N0RK5/ I NUXIH  10) 


THIS  SUBROUTINE  FORMS  MIC  13  AND  C MATRICES  BASED  UPON  THE  DESIGN 


I F ( IPOS.EO.O)  GO  TO  13 
I F ( IPOS.EQ.2)  GO  TO  d 
DO  0 11*1, IUNKR 
I * I NDXK( I I ) 

RRC=K( I ) *RC(  I ) 

IFURC.GE.u.O)  RRC=6 . 0 
P( I )*EXP(RRC) 

PC(l)=i.O 

I F ( ITRAN.EO.O)  P( I ) *R( I ) 

I F(  I TRAii.E  1.0)  PC  ( I ) =RC  ( I ) 

CON  I I NUE 
GO  TO  13 

00  12  II *1, IUNKR 
I =INOXR( I I ) 

IF(P( I ) .LE.l.E-3)  GO  TO  10 

1 F ( .1C ( I ) . LE  .1 . E-8 ) RC(  I )-1.0 
R(  I ) = ALOG(  l'(  I ) )/RC(  I ) 

GO  TO  11 
RCI )-  -10.4 
CONTINUE 

Id  I TRAN.EQ.O)  R ( I ) = P(  I ) 

I F( I TRAN. EG. 0)  RC(I)=PC(I) 

CONTINUE 
CONTINUE 
UCOM*14 . Oil) 

0011*1, NA 
D01J=1, NA 
0(l,J)=0.0 
CONTINUE 
0021=1, NA 
002 J=1,MC 
C( I , J ) ” 0 . 0 
CONT I NUE 

20  = -1.0*(()U.  4 6 50  3 7 *UCOM )*P(2)*PC(2)*(-22. 2 692  70  *UC0M)  » 

1P(3)* PC (5)*(lll. 627245 *UCO!l) ) 

U(l,l)=-1.0*((l. 10  2744  *UC0M) >P( 1)*PC(1)*(0.075S73*UCOM)  + 
1P(2)*PC(2)*(0. 0441O3  *UCOM) ♦P(3)*PC(3)*(U. 074  750*UCOM) *P(4)*PC(4)* 
2 ( 0 . do  1020  *UCOM) ) 

6(1, 1)=-21.0340 
a(l,3)=0.QWQ 
6(1, 4) =5. 3721 

B(2,2)=-1.0»((3.6‘J323*UCOM)+P(1)*PC(1)»(0 . 037743  *UC0M) ♦P(2)*PC(2) 
1*( 4 . 7o4  3o J*UCOM)  + P( 3)*PC(3)*(6.01 750*UCOM) ♦P(4)»PC(4)* 

2(0. 030u006*UCOM) ) 

B(2, 3)=ZQ+200. 3d*UC0M 

B(5,Z)=(iJ2. 156*UC0M) *P( 2 ) *PC( 2 ) * ( 25 . 09S000*UCOM)+P( 3) *PC( 3) * 


34  - 


r 


21u2d 
2 lu  3 J 
2 1 u >»  J 
21udJ 
2 lu  j J 
2 lo  / J 
2 looO 
21uJd 
21/JJ 
2 1 / id 
2 1 1 2d 
21/ jj 
2i/sd 
2 1 / j d C 
21/od  C 
2 1 / / d C 
2 1 / d J 
21/JJ 

2 lo  d d 
2 lo  Id 

2 1 o 2 d 
2 JLo  3 d 

2 1 o •*  J 
2 lo  J J 

2 1 o u d 
21o/d 

2 lo  d d C 
2 Id  Jd  C 
21 J JJ  C 
21  Jld 
2 1 J 2 d 
2l  J id 
2 1 J 4 J 
2 1 J j J 

2 1 Jo J 3 
21J7J  C 

2 1 J d J C 
21JJJ  C 

22  J J J 
22dld 
22  J 2d 

2 2 J 3 0 ** 
22  J4d  C 
2 2 J J J C 
2 2 J o J C 
2.10 1 J C 
22000 
22000 
22000 
22000 
22020  C 
2213d  C 
22140  C 

221  Jd 
221o  J 
2 i 1/0 

221od 

221JJ 
222 Jj  j 

2221J  C 
222i0  C 
2223d  C 
2 2 24  J 

222  jd 
2 2 2 od 
22270  o 
222dd  C 
222  JJ 
22303 


f*s?AG1; 

IKOV 


1( -35 . 5 179*UCOM> 

!i(3, 3)  = -JJl. 436809  *UC0M»P(  2)*PC(2)*(-10S.4  5 114  3 *UC0lt)  * 
1P(3)*PC(3)* ( -1596 . 29dd  3 7*UC0M) 

U( 3, 4 )=- 722 . 2 
B<4,3)=1.0 
d( j, 2 )*1 . d 
B( j,4)=-UC0M 

C(1,1)=-P(2)*PC(2)*(0.55'J1I)6)*(  UC0M*  *2 ) 

C( 1, 2 ) =-P( 3) *PC( 3 ) * ( U . 250559) *( UC0M* * 2 ) 

C ( 2,  1 ) =-P(  2 ) *PC(  2 ) * ( 3 . 5 J lU5u ) *(UC0il**  2 ) 
C(2,2)=-P(3)*PC(3)*(2.5d5594) *( UC0M**2 ) 

C(5, i)=P(2)*PC(2)*(27.223202)*(UC0M*»2) 
C(J,2)=P(3)*PC(3)*(-3/.2d3137)*  (UCOIt*  *2 ) 


FORM 


E MATRIX 


£(l,l)=d.d 
E(2,l)=d.O 
E( 1, 2 ) =P( 5 )*PC( 5 ) 

E(i,3)=-P(o)*PC(6) 

E ( 1 , 4 ) =- P( lo)*PC(16)*P(d)*PC(8)-P( j)*PC(5)*UC0M 
E(1,5)=P(15)*PC(15)*P(7)*PC(7) 

E(2,2)=-P(9)*PC(9) 

E(2,3)=P(10)*PC(1J) 

E(2,4)=P(14)*PC(14)*P( 12)*PC(12)+P(9)*PC( J) *UC0M 
E(2,5)=-P(15)*PC(13)*P(11)*PC(11) 


MULTIPLY  C AMD  E MATRIX 


00  31=1,  !4A 
D05J=1,NA 
CE( I , J)=J.U 
00  3K=1,MC 

CL( I , J ) =CE ( I , J ) ♦ C(I,K)*£(K,J) 
CONTINUE 


FORM  OVERALL  B MATRIX 

00  41=1, NA 
00  4j=l,NA 

B ( I , J ) = d ( I , J ) ♦ C E ( I , J ) 
CONTINUE 

FORM  0 MATRIX 


DU,  1)=P<S)*PC(3) 
U(i,2)=-P(7)»PC(7) 
D(2,1)=-P(12)*PC(12) 
D(2,2)=P(11)*PC(11) 

FORM  NE.J  CONTROL  MATRIX 

DO  51=1, NA 
00  5 J = 1 , MC 
CE ( I , J ) = J . 0 
DO  5 K = 1 , MC 

CL( I , J ) =CE( I , J ) ♦ C(I,K)«D(K,J) 
CONTINUE 

FORM  OVERALL  C MATRIX 

00  ul = 1 , UA 
00  o J = 1 , MC 
C( I , J) =CE( 1,0) 

CONTINUE 

RETURN 

END 
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